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Possible  Estimation  Methodologies  for  Electromagnetic 
Field  Distributions  in  Complex  Environments 

M.  Kanda,  J.  Randa,  and  N.  S.  Nahman 

Electromagnetic  Fields  Division 

National  Bureau  of  Standards 

Boulder,  Colorado  80303 


The  problem  of  measuring  and  characterizing  complicated  multiple- 
source,  multiple-frequency  electromagnetic  environments  is  becoming 
more  important  and  more  difficult  as  electrical  devices  prolifer- 
ate. This  paper  outlines  three  general  approaches  to  the  problem 
which  are  currently  under  investigation  at  the  National  Bureau  of 
Standards.  The  three  approaches  are:  1)  a  statistical  treatment  of 
the  spatial  distribution  of  electromagnetic  field  intensities,  2)  a 
numerical  computation  using  a  finite-element  (or  lattice)  form  of 
the  electromagnetic  action  functional,  and  3)  use  of  a  directional 
probe  to  scan  a  volume.  All  three  methods  are  still  in  the  develop- 
ment stage,  but  each  appears  promising. 

Key  words:  action,  directional  scanning,  environment  characteriza- 
tion, field  levels,  finite  element,  hazard  assessment,  multiple 
source,  statistical  approach. 

1.  Introduction 

There  has  recently  been  a  rapid  increase  in  higher  powered,  multi- 
frequency,  electromagnetic  (EM)  radiation  sources  which  complicate  greatly  the 
environment  in  which  modern  electronic  equipment,  both  military  and  civilian, 
must  operate.  At  the  same  time,  there  has  been  an  even  more  dramatic  increase 
in  the  quantity  of  electronic  equipment,  such  as  minicomputers  and  micro- 
processors, composed  of  semiconductor  devices,  and  hence  more  sensitive  to 
interference.  Today  there  is  a  complex  matrix  of  electronic  equipment  trying 
to  operate  in  a  complex  EM  environment.  Therefore,  the  estimation  of  the  max- 
imum electromagnetic  (EM)  field  strength  is  becoming  very  important  in  many 
interference  problems.  For  example,  electromagnetic  waves  penetrate  into 
buildings  which  house  sensitive  electronic  and  ordnance  items.  There  has  been 
great  interest  in  determining  the  electromagnetic  environment  inside  build- 
ings. 

Heretofore,  it  has  been  customary  to  make  several  spot  measurements  of 
the  electromagnetic  field  strengths  produced  by  existing  transmitters.  How- 
ever, it  is  obviously  impossible  to  completely  determine  the  field  existing  in 


an  enclosure  using  a  manageable  number  of  spot  measurements.  In  order  to 
rigorously  determine  the  EM  field  distributions  in  a  volume,  systematic  mea- 
surements of  electric  field  amplitude  and  phases  have  to  be  made  at  points  at 
most  half  a  wavelength  apart  (at  the  highest  frequency  present). 

There  is  a  need  for  a  general  method  of  extracting  the  maximum  amount  of 
useful  information  about  the  EM  field  distribution  within  an  enclosure  from 
the  minimum  amount  of  measurement  effort.  In  this  paper  we  discuss  three  pos- 
sible approaches  to  this  problem  which  are  currently  being  explored  and  devel- 
oped at  the  National  Bureau  of  Standards  (NBS) .  We  first  examine  a  statisti- 
cal approach  to  the  estimation  of  EM  field  distributions.  In  a  previous  paper 
by  the  first  author  [1],  five  types  of  time  and  amplitude  statistics  were  used 
in  order  to  unravel  the  complexities  involved  in  an  EM  environment.  Here  we 
shall  discuss  the  statistical  distribution  of  scattered  EM  fields  and  illus- 
trate its  usefulness  using  experimental  and  simulated  examples.  The  second 
method  is  a  finite-element-action  calculation  using  a  small  number  of  measure- 
ment points.  The  true  field  cannot  be  determined,  but  it  may  be  possible  to 
reconstruct  the  smoothest  configuration  or  the  "most  probable"  one.  The  third 
and  final  approach  discussed  is  the  most  conventional,  using  directional 
scanning  to  bound  or  approximately  determine  fields  within  a  given  volume. 

2.  Statistical  Approach 
2.1  Rayleigh  Distribution 

When  the  configuration  of  the  scatterers  is  random  and  sufficiently  dis- 
persed to  give  a  wide  (and  hence,  equivalently  uniform)  phase  distribution, 
the  field  configuration,  E,  typically  caused  by  the  reflectivity  of  various 
multipath  types  of  scatterers,  both  stationary  and  moving,  is  given  by 

E  =  I      A.ej*i,  (2.1) 

i=l   1 

where  A^  and  <j>^  are  random  amplitudes  and  phases  of  scattered  fields. 

Resolving  E  into  its  real  and  imaginary  components,  Re  E  and  Im  E,  from  the 
central  limit  theorem,  it  follows  that  Re  E  and  Im  E  are  normally  distributed, 
as  long  as  N  is  large  (>  10).  In  particular,  when  the  A.,-  and  ^  are  uncor- 


related  and  hence  independent  random  variables,  it  can  be  shown  that  Re  E  and 

Im  E  have  the  same  variance  with  zero  mean  values.  The  amplitude  distribution 

r       r 
pF(r)  is  then  found  to  be  the  Rayleigh  distribution  [2]  pr(r)  =  — ~  exp  (-  -r— 5) 
D  E  a  do 

for  r>0. 

The  mean  <E>,  mean  square  <E^>  and  standard  deviation  (s.d.)  are  respec- 
tively given  as 


<E>  =  ]  rpF(r)  dr  =  /  jr  a  (2.2) 

0    L         2 

<E2>  =  2a2  (2.3) 


andv 


s.d.  =  /<E2>  -  <E>2  =  /2  -  -  o.  (2.4) 


The  cumulative  probability  distribution  is  defined  as  the  probability  that  the 
random  variable  E  is  equal  or  less  than  the  value  E0  ,  and  can  be  obtained  by 
integrating  the  Rayleigh  probability  density  function,  i.e., 

Eo  c  2 

2  0 

P(E  <  E  )  -  J  4  exP  ("  Tl)   dr  "  1  "  exP  C"  TT)-        (2-5) 
0    q  a  co  co 

Hence,   its  complement  is 

Eo2 
p   (E    >  E0)   =  exp(-  y^i)  (2.6) 

The  median  of  E,   Em,  is,  therefore, 

E  2 

P    (E  <Em)  =  1   -  exp    (- JJJ2)  =  0.5  (2.7) 

Em   s  1.18   a.  (2.8) 


or 


The  average  crossing  rate  (ACR)  presents  the  average  number  of  times  the 
EM  field  strength  crosses  various  levels  and  is  usually  given  as  positive 
crossing  per  second  versus  EM  field  strength.  The  average  level  crossing  rate 
at  a  given  level  EQ  is  [3] 

00 
n(E  =  E0)  =  /  p  pE(E0;r)  dr  (2.9) 

0 


where  pF(Eo,r)     is  the  joint  probability  distribution  of  E,  and  p  =  -rp    gives 
the  slope  of  the  signal   envelope  E.     For  the  case  of  the  Rayleigh  distribu- 
tion, the  level   crossing  rate  becomes 

E  F   2 

v         °  ° 

n   (E  =  E0)   =  — -    -^    exp(-  -^)  (2.10) 

where*  *^n 


a2    **\>\r  2ai. 


v2  =  <E2>  .  (2.11) 

The  average  field  strength  can  be  determined  by  averaging  the  random 
trial   samples  of  instantaneous  field  strength  measurements   recorded  over  a 
certain  path  length  in  distance.     The  variation  in  the  average  field  strength 
is  typically  caused  by  the  relatively  small   scale  variations  along  the  propa- 
gation path. 

There  are  situations  where  a  signal   is  received  in  two  different  ways, 
one  via  many  widely  spaced  scatterers  or  reflectors,  the  other  by  a  quite  dif- 
ferent mechanism,  e.g.,  directly  from  a  transmitter.     The  analysis  of  the  am- 
plitude distribution  of  the  resultant  signal   can  be  achieved  as  follows.     Con- 
sider the  sum 

•a      N 
E  =  A0eJ<l>+     I      A.eJ<|,i.  (2.12) 

i=l 

The  first  term  represents  the  direct-path  propagation  through  the  walls  of 
houses  and  other  objects.  The  second  term  represents  the  reflected  and  scat- 
tered waves  due  to  randomly  oriented  objects  such  as  houses,  buildings,  walls, 
and  overhead  wires,  and  is  therefore,  Rayleigh  distributed.  Since  the  direct- 
path  propagation  is  attenuated  through  each  penetration  of  a  wall,  etc.,  the 
amplitude  of  this  direct  wave  will  be 

M 
Ao  ttexp  [-  I     6.  d.]  (2.13) 

i=l  ]  ^ 

where  d^  is  the  thickness  and  5j  is  the  attenuation  constant  of  each  of  the 
walls.  When  M  is  large,  the  central  limit  theorem  indicated  that  the  sum  of 
the  random  term  6.,-d.j  will  be  normally  distributed.  Hence,  the  amplitude  dis- 
tribution of  A0  is  lognormal,  i.e., 


Pn  -(r)  = 


1 


r  (Inn  -  y)^ 


(2.14) 


0     r  a  /2tt 
where  y  is  the  mean  and  a2  is  the  variance  of  the  exponent  in  eq  (2.13). 

The  probability  density  function  of  the  total  field  given  by  the  sum  of  a 
random  plus  a  Rayleigh  phasor  is  [2] 


P(r)  = 


2r 


a  a  /2ir  ° 


/  A"0  exp[- 


(*n  An  -  y)2   r2  +  A  2 


2r  A, 


2  a2 


]  -In  (" 


■)  dA 


(2.15) 


where  a  is  the  mean-square  value  of  the  scattered  components. 

For  the  cases  of  large  and  small  values  of  E, 

1  2  - 

exp    [-   ( jm  r  -    y)   /2a2]     for  r  >  /2a 


PE(r)   =  < 


ra/2ir 


exp    (-   r2/2o2) 


for  r   <  /2a 


where   y  is  the  mean  and   a2  is  the  variance  of  E. 


(2.16) 
(2.17) 


Equations  (2.16,  2.17)  indicate  that  E  is  lognormal  for  large  values  and 
Rayleigh  distributed  for  small  values. 

The  Rayleigh  distribution  has  one  parameter  a  to  be  determined  which  fol- 
lows  simply  from  the  estimated  mean  a  given  as 


1 


N 


and  therefore 


P=  N  i=Zi  Ei' 


(2.18) 


a   =  /  —  y. 


(2.19) 


2.2.  Examples 

A  number  of  electromagnetic  field  measurements  have  been  made  on  the 
propagation  of  radio  waves  at  VHF  and  higher  frequencies  in  the  presence  of 
buildings,  trees  and  other  obstacles  [1,4,5,6].  As  an  example,  the  cumulative 
amplitude  probability  distribution  and  the  average  crossing  rate  of  the  EM 
noise  measured  in  a  coal  mine  [1]  are  shown  in  figures  2-1  and  2-2, 


respectively.  Figure  2-1  indicates  that,  when  the  field  strength  is  low 
(approximately  90%  of  the  total  sample),  the  Rayleigh  distribution  is  a  good 
approximation  of  the  EM  noise  in  a  mine.  On  the  other  hand,  when  the  EM  field 
strength  is  high  (less  than  10%  of  the  total  sample),  EM  noise  data  show  a 
departure  from  the  Rayleigh  distribution.  The  statistical  properties  of  the 
peak  values  of  EM  field  may  be  close  to  the  log-normal  distribution  as  indi- 
cated in  eq  (2.16) . 

Measurements  of  attenuation  due  to  buildings,  trees,  and  other  obstacles 
at  long  propagation  distances  (e.g.,  ranging  from  100  to  500  m  for  800  MHz 
propagation),  have  also  been  performed  [6],  As  an  example,  figure  2-3  shows 
the  cumulative  amplitude  probability  distribution  of  the  attenuation  data.  A 
straight  line  on  the  coordinates  in  figure  2-3  represents  a  log-normal  distri- 
bution of  signal  levels.  The  measured  distribution  is  approximated  by  a  log- 
normal  distribution  which  is  discussed  in  eq  (2.16). 

To  simulate  these  experimental  results,  numerical  studies  of  the  cumula- 
tive amplitude  probability  distribution  for  a  two-dimensional  cavity  model 
with  random  noise  sources  have  been  made  [7],  Figure  2-4  shows  the  cumulative 
amplitude  probability  distribution  for  EM  fields  measured  inside  the  rectangu- 
lar cavity.  This  figure  is  obtained  from  100  data  sets  which  are  generated  by 
using  random  noise  source  positions  and  random  wavelengths.  The  results  fol- 
low a  lognormal  distribution  except  for  the  low  probability  region  where  most 
electromagnetic  interference  data  show  a  departure  from  a  lognormal  distribu- 
tion. The  deviation  from  a  lognormal  distribution  in  the  low  probability 
region  indicates  an  accurate  estimation  of  the  maximum  field  strength  for  the 
measured  data  may  be  very  difficult. 

2.3   Other  Relevant  Statistics 

In  its  most  general  form  the  characterization  of  a  stationary  electro- 
magnetic environment  is  based  on  observations  of  a  three-dimensional,  random 
variable.  As  such,  the  electromagnetic  field  is  represented  by  the  electric 
and  magnetic  field  vectors  E   and  TT,  respectively,  where  both  of  these  vectors 
may  be  stationary  random  functions  of  their  spatial  coordinates. 

If  the  field  vectors  are  not  stationary,  the  probability  distributions  of 


the  electromagnetic  field  vectors  are  dependent  upon  time.  The  process 
required  for  the  environment  to  change  from  one  stationary  state  to  another  is 
defined  as  the  transient  state  and  is  not  considered  here. 

For  stationary  electromagnetic  excitation  the  applied  fields  or  sources 
whose  collective  effects  establish  the  EM  environment  within  a  volume  are 
defined  here  as  steady-state  fields.  Thus,  the  random  nature  of  the  field  is 
established  by  the  random  spatial  distribution  of  the  sources,  and  the  random 
boundaries  of  objects  within  the  volume  or  the  randomness  of  the  boundary 
itself. 

There  are  various  statistical  approaches  which  may  be  useful  for  charac- 
terizing this  class  of  electromagnetic  environments;  two  approaches  are  pre- 
viously discussed;  i.e.,  cumulative  probability  distribution  and  averaging 
crossing  rate. 

Three  other  approaches  are  mentioned  here  as  having  possible  merit. 

A.  Interpulse  Spacing  Distribution 

The  interpulse  spacing  distributions  give  the  probability  distribution 
for  the  spacing  between  successive  pulses  in  the  received  noise  process. 
These  distributions  are,  of  course,  functions  of  the  noise  amplitude  level. 

B.  Pulse  Duration  Distribution 

The  pulse  duration  distributions  give  the  probability  distribution  for 
the  pulse  widths  and  are  given  in  terms  of  the  percentage  of  pulses  which 
exceed  various  widths  in  seconds. 

C.  Two  Sample  Variance  Analysis 

It  is  essential  to  know  how  much  data  to  gather  when  dealing  with  statis- 
tical quantities.  Therefore,  in  any  measurement  of  a  statistical  phenomenon 
the  minimum  length  of  time  over  which  the  phenomenon  is  observed  should  be 
determined.  Two  sample  variance  analysis  can  be  used  to  accomplish  this 
determination.  The  basic  idea  to  be  discussed  briefly  below  has  been  imple- 
mented often  in  the  discussion  of  frequency  stability. 


A  record  of  the  phenomenon  under  consideration,  y(t),  is  divided  into  a 
number  of  equal  time  segments  of  length  t,  and  the  average  value  of  y(t),  y^, 
of  each  segment  is  calculated  by 


VT 


=  ±   j   y(t)  dt,  (2.20) 


yk 

t 

where  y^  is  the  kth  segment  average  starting  at  time  t^.  Next,  the  sample 
variance  (sample  size  two)  o2  (2,t),  of  successive  averages  is  calculated. 
That  is 


where 


0kM2,,)=    "?   (yn->"k)2=|(Vk+r^)2-  <2-21> 

n=k 

1   k+1 
yk  4   e  yn  (2.22) 

n=k 

is  the  average  of  the  two  successive  segment  averages  y^  and  y^.  The  two 
sample  variance,  a2  (2,t),  for  this  special  case  (sample  size  two)  is  then 
defined  to  be 

a2  (2,t)  =   <a2  (2,t)>,  (2.23) 

where  the  brackets  represent  the  average  of  a2  (2,t)  over  all  pairs  of  suc- 
cessive yk  constructed  from  y(t).  The  preceding  calculation  is  repeated  for 
various  values  of  averaging  period,  t.  For  a  given  maximum  allowable  devia- 
tion in  y(t)  the  minimum  averaging  time  can  then  be  determined. 

Many  examples  of  these  time  and  amplitude  statistics  for  the  time  depen- 
dent, EM  noise  are  given  elsewhere  [1], 

3.     Finite-Element  Action  Approach 
3.1      Introduction  and  Motivation 

In  this  section  we  describe  an  approach  which  in  effect  attempts  to  solve 
Maxwell's  equations  within  the  volume  of  interest.     The  type  of  problem  in 
which  we  are  interested,  however,  differs  fundamentally  from  those  problems 
for  which  Maxwell's  equations  are  usually  solved,  and  it   resists  the  use  of 
standard  methods.     The  most  fundamental   difference  is  that  we  do  not  know  the 


sources.  In  addition  there  are  the  complications  that  the  geometry  and  the 
time  dependence  are  not  simple.  The  information  available  includes  boundary 
conditions  at  conducting  walls  or  dielectric  interfaces  within  the  volume,  and 
the  measured  values  of  the  fields  at  some  number  of  measurement  points.  We 
are  free  to  specify  the  number  and  location  of  the  measurement  points,  but  we 
obviously  want  the  number  to  be  small — or  at  least  "reasonable."  The  general 
qualitative  idea  is  that  given  the  geometry  and  boundary  conditions,  knowledge 
of  the  fields  at  a  few  points  should  enable  one  to  extract  some  information 
about  global  properties  such  as  average  or  maximum  electromagnetic  energy  den- 
sity. 

Without  knowledge  of  the  external  sources  radiating  into  the  volume  of 
interest,  and  with  only  a  few  measurement  points,  there  will  not  in  general  be 
a  unique  solution  to  Maxwell's  equations.  Consequently,  it  will  be  impossible 
to  actually  determine  the  field  everywhere.  The  best  one  can  hope  for  is  that 
the  approximate  values  and  positions  of  the  maximum  field  intensities  will  be 
the  same  for  all  allowed  solutions.  Failing  that,  one  wants  a  method  which 
consistently  finds  the  one  solution  of  all  those  possible  which  has  some  prop- 
erty of  interest,  such  as  the  smallest  maximum  power  density. 

An  approach  which  appears  to  hold  some  promise  of  fulfilling  these  rather 

demanding  (perhaps  impossible)  requirements  is  based  on  a  finite-element 

treatment  of  Hamilton's  principle  applied  to  electrodynamics.  For  a  classical 

field  theory,  Hamilton's  principle  states  that  if  one  considers  the  quantity 

t2       t2 

S[*a]  =  J      dt  L  =  /      dt     |d3x    t(*a(x,t)),  (3.1) 

ti  tx 

where  £.(\|>  )   is  the  Lagrangian  density  of  the  system  depending  on  the  indepen- 
dent fields   ij>a(x,t),  then  the  values  assumed  by  the  fields  for  the  correct 
physical   solution  are  such  that  S  is  stationary  with   respect  to  small    varia- 
tions of  the  fields,   V*(x,t).     The  functional    S[^a]  is  called  the  action,  and 
for  electrodynamics  the  stationarity  requirement  yields  Maxwell's  equations 
[8],  as  we  shall   outline  in  the  next  subsection.     It   is  generally  assumed  that 
the  stationary  point  of  the  action  is  in  fact  a  minimum,  and  general   condi- 
tions are  known  for  which  this  is  so  [9],     Unfortunately,  these  conditions  are 
generally  not  met  in  the  cases  of  interest  to  us,  and  so  finding  the 


stationary  point  of  S  is  not  necessarily  as  "simple"  as  minimizing  it.  This 
will  be  discussed  further  in  a  later  subsection. 

Finite-element  calculations  minimizing  the  action  have  been  used  before 
in  electromagnetics,  particularly  for  waveguide  and  cavity  problems  [10,11], 
and  the  connection  to  other  variational  calculations  has  been  pointed  out 
[11].  Previous  approaches,  however,  have  obtained  a  set  of  linear  equations 
by  using  finite  elements  and  setting  the  appropriate  derivatives  equal  to 
zero.  The  system  was  then  solved  by  numerically  inverting  the  matrix.  Such 
tatics  will  not  work  here  because  there  is  not  a  unique  solution,  and  hence 
the  matrix  would  be  (very)  singular.  Our  preference  is  to  work  directly  with 
the  action  and  search  numerically  for  its  stationary  points,  rather  than  work- 
ing with  its  derivatives.  This  approach  is  probably  closer  in  spirit  to 
methods  currently  popular  in  quantum  field  theory  [12]  than  to  the  traditional 
differential  equations  approach. 

In  the  finite-element  approach  we  replace  the  continuous  variables  x  and 
t  by  a  four-dimensional  grid  of  points  on  which  the  fields  are  defined.  The 
integral  in  the  action  then  becomes  a  sum,  and  the  field  values  at  each 
(unmeasured)  point  are  varied  until  a  stationary  point  of  the  action  is 
found.  The  field  values  at  measurement  points  are  set  equal  to  their  measured 
values  and  not  allowed  to  vary.  When  a  stationary  point  is  found,  one  has  a 
solution  of  Maxwell's  equations  for  which  the  fields  take  on  their  measured 
values  at  measurement  points.  Which  of  the  many  possible  such  solutions  one 
finds  will  depend  on  the  initial  field  configuration  (before  we  vary  the 
fields)  and  on  the  method  used  to  find  the  stationary  point.  In  the  present 
work  we  are  using  a  smooth  starting  configuration  and  a  gradient  "minimiza- 
tion" procedure,  which  should  find  the  solution  which  has  the  smoothest  field 
configuration. 

The  remainder  of  this  section  is  devoted  to  a  more  detailed  account  of 
the  progress  made  in  formulating  this  approach  and  the  (foreseeable)  remaining 
problems.  We  shall  first  present  the  general  framework,  then  digress  on  the 
question  of  whether  the  stationary  point  is  a  minimum.  Next  a  simple  example 
is  given,  and  finally  we  discuss  remaining  obstacles  and  the  direction  of 
future  work. 
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3.2  General  Formulation 

The  action  for  electromagnetism  is  given  by  [8] 

S[A,*]  =  ]   dt  ]  d3x  {i  (V<|>(x,t)  +.|j-A(x,t))  •  e(x,t)  .  (v*(x,t) 
1 1 
+  -^-A(x,t))  -  (V  x  A(x,t))  •  u_1(x,t)  •  (V  x  A(x,t)J] 

+  A(x,t)  •  J(x,t)  -  p(x,t)  ♦(x,t)}  ,  (3.2) 

where  we  have  assumed  all   media  are  lossless,  and   e  and   y  are  the  permittivity 
and  permeability  tensors.     A  and   <|>  are  the  usual   vector  and  scalar  potentials, 

__  _3_  __  __ 

E(x,t)    =  -V<fr  -  —  A        ,       B(x,t)    =   V  x  A      .  (3.3) 

ot 

The  time  integral   in  eq   (3.2)   is  typically  taken  from  -»  to  +  °°,  and  similarly 
the  volume  integral   usually  extends  over  all   space--al though  one  can  also  con- 
sider a  finite  volume  or  time  range,  subject  to  conditions  we  shall   note 
below.     J  and  p  are  the  current  and  charge  densities   respectively.     We  have 
written  the  action  in  terms  of  the  potentials   rather  than  E  and  H,  because  in 
varying  the  action  we  want  to  take  variations  only  with   respect  to  independent 
variables,  whereas  E  and  H  are  related  through  eq   (3.3).     In  fact,  eq   (3.3) 
implies  the  two  homogeneous  Maxwell's  equations, 

V-B  =  0       ,       VxE+  —  B  =  0.  (3.4) 

3t 

Although  there  are  four  independent  potential  functions  (A,<|>)  appearing  in  eq 
(3.2),  only  three  are  independent,  since  the  gauge  choice  will  provide  an 
extra  relation. 

If  one  considers  variation  of  the  action  due  to  a  small  variation  of 
<t>(x,t)  (6<j)(x,t))  one  finds  (assuming  the  permittivity  tensor  is  symmetric) 
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t2  t2 

6  S  =   -   /       dt     /     d2x  n    •  D   6<j>(x,t)   +     /       dt   /  d3x   (V   •  D  -   p)    6<j>, 

*i  s  li  (3.5) 

5  -   e    •   (-V*  --JJ-A), 

where  the  surface  integral   is  over  the  surface  enclosing  the  entire  volume. 
If  we  require  that   <J>  not  vary  on  the  boundary  surface,   but  otherwise  let 
6cj>(x,t )   be  arbitrary,  then  the  first  term  in  eq   (3.5)   vanishes,  and  requiring 
.that  the  action  be  stationary  yields 

6  S   =   /       dt      j  d3x    (V   •  D   -   p)    6<j>  =   0 

*l  (3.6) 

*  V   •  D(x,t)    -   p(x,t)   =  0   . 

In  a  similar  manner,   requiring  that  S  be  stationary  with  respect  to  small 
variations  of  A(x,t)   leads  to 

6S=0     +    —  6-VxH  +  J  =  0,  (3.7) 

A  dt 

where  A(x,t)  must  be  held  fixed  not  only  on  the  boundary  surface  at  all  times 
but  also  throughout  the  volume  at  times  tj  and  t2»  What  this  means  in  prac- 
tice is  that  in  order  to  use  this  in  a  calculation,  one  must  specify  the 
fields  on  the  boundary  surface  at  all  times  and  throughout  the  volume  at  the 
initial  and  final  times.  In  any  real  application,  such  a  superabundance  of 
information  will  not  be  available  for  the  volume  of  interest.  In  order  to 
exploit  the  stationarity  of  the  action  we  must  expand  the  volume  considered 
beyond  just  the  region  of  interest,  out  to  distances  where  the  fields  can  be 
assumed  to  be  negligible.  The  same  applies  to  the  time;  tj  is  chosen  before 
the  fields  are  turned  on  and  t2  after  they  are  turned  off. 

Having  extended  the  volume  under  consideration  to  virtually  all  of  space- 
time,  we  must  restrict  the  problem  to  manageable  size.  We  envision  the  divi- 
sion depicted  in  figure  3-1.  The  volume  marked  Vj  is  the  region  of  interest; 
it  is  assumed  to  be  free  of  primary  sources,  but  it  may  contain  conductors 
with  induced  currents.  Volume  VB  is  a  buffer  zone  between  Vj  and  Vs,  which  is 
the  rest  of  space,  wherein  are  located  any  primary  sources.  The  general  idea 
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is  to  divide  the  action  into  one  piece  from  the  integral   over  volumes  Vj  and 
Vg  and  another  piece  from  V<;.     Because  V^  contains  unknown  sources,  we  do  not 
try  to  determine  the  fields  there,  which  means  that  we  also  will   not  know  the 
fields  on  the  surface  between  Vs  and  Vg.     The  fields  on  that  surface  will   be 
allowed  to  vary  or  will   be  fixed  by  a  reasonable  guess.     Obviously,  near  the 
boundary  between  V<j  and  Vg  the  solution  obtained  will   be  very  sensitive  to  the 
choice  for  the  fields  on  that  boundary,  and  therefore  it  will    not  be  reli- 
able.    As  one  moves  away  from  the  outer  boundary  of  Vg  the  values  of  the 
fields  should  be  influenced  more  by  the  measurement  points  and  less  by  the 
values  on  the  surface  between  V$  and  Vg.     For  points  far  enough  away  from  that 
surface--i  .e.  within  Vj — it  should  be  possible  to  obtain  reliable  solutions, 
given  enough  measurement  points.     The  hope  is  that   "far  enough  away"  and 
"enough  points"  are  not  so  large  as  to  render  the  method  impractical   for  most 
applications.     The  positions  of  the  measurement  points  will   clearly  affect  the 
size  required  for  Vg;   it  may  well   be  advantageous  to  make  a  few  measurements 
on  the  perimeter  of  Vg.     It  would  probably  also  be  advisable  to  choose  the 
boundary  of  Vg  to  coincide  with  conducting  walls  when  it  is  feasible,   in  order 
to  constrain  the  fields  on  the  boundary  as  much  as  possible. 

The  quantity  we  shall   consider  then  is  a  reduced  action,  S,  which   is 
defined  as  in  eq   (3.2)   but  with  the  integral    restricted  to  V.  +  VR.     (Note 
that  in  the  general   multiple-frequency  case  the  various   volumes  are  four- 
dimensional   space-time  volumes.)     For  definiteness  and  simplicity  we  choose  a 
gauge  which  will   be  used  in  the  remainder  of  the  section.     The  choice  is 
<|>(x,t)   =  0,  so  that  E  =  -3A/3t .  The  next  step  is  to  discretize  the  expression 
for  the  reduced  action,  converting  the  volume  integral    into  a  summation  which 
approximates  it.     There  are  any  number  of  ways  to  do  so;  we  are  not  interested 
in  their  relative  (dis)advantages  at  this  time.     The  discretized  reduced 
action  will   have  the  general   form 

S  =     /  dt  I  AV         (i(^A(t))  •  I  v(!rA(t)) 

a.B.YeVjfcVg  (3>8) 

-    (V  x  A(t)l  •   y"1      •    (V  x  A(t))        }   , 

where  a  number  of  points  require  explanatory  comment.  The  discrete  indices 
a,3,Y  label  the  spatial  points,  on  which  are  centered  the  volume  elements 


13 


AV  „  .   We  have  left  the  time  variable  continuous  for  now,  anticipating  the 
<*3y  r    3 

single-frequency  example  below,  Section  (3.3),  but  for  the  general  case  time 

too  would  be  made  discrete.  The  quantities  (8A/8t)    and  (V  x  A)    will  in 

a3Y  a&Y 

general  depend  on  the  values  of  the  field  A  at  a  number  of  points  on  (or 

within)  the  surface  bounding  AV  _  .  A  concrete  realization  of  this  discretiza- 
*  3    cx3Y 

tion  of  S  will  be  given  in  the  example  below.  Imposition  of  the  constraints 

required  by  the  known  values  of  E  and/or  H  at  measurement  points  can  be  rather 

complicated,  but  a  simple  case  illustrates  the  idea.  If  only  one  frequency  is 

present  then  for  our  gauge  choice  E    «  A   ,  and  a  measurement  of  the  elec- 

a3Y     a3Y 
trie  field  at  a  point  fixes  A    at  that  point.  A  similar  comment  applies  to 

a3Y 
boundary  conditions  at  perfectly  conducting  walls;  for  the  single-frequency 

case  they  can  be  imposed  with  relative  ease. 

The  calculation  then  proceeds  as  follows.  A  grid  is  defined  within 

Vg^Vj,  and  A    is  fixed  at  measurement  points  and  appropriate  components  are 

set  equal  to  zero  at  conducting  walls.  All  other  A   's  are  varied,  and  one 

a3Y 
searches  numerically  for  a  stationary  point  of  S  in  eq  (3.8).  When  (if)  a 

stationary  point  is  found,  then  that  set  of  A    constitutes  an  approximate 

cx3y 
solution  to  Maxwell's  equations  within  Vg  ®  Vj  which  is  consistent  with  mea- 
sured field  values  and  boundary  conditions.  The  entire  procedure,  including 
the  numerical  search  method,  is  greatly  clarified  by  an  example,  to  which  the 
next  subsection  is  devoted.  The  question  of  whether  the  stationary  point  is  a 
maximum,  minimum,  or  neither  is  addressed  in  the  Appendix. 

3.3  Simple  Example 

Having  presented  the  general  ideas  of  this  approach,  we  now  attempt  to 
implement  it  in  a  simple  example — a  rectangular  waveguide  with  perfectly  con- 
ducting walls.  This  is  obviously  not  supposed  to  be  a  practical  application, 
and  many  of  the  difficulties  and  nuances  of  the  general  case  are  absent.  It 
is  a  practice  problem  to  demonstrate  the  idea  and  to  provide  a  basis  on  which 
to  build  toward  solution  of  real  problems.  The  rectangular  waveguide  is 
chosen  because  it  has  simple  known  solutions,  because  the  boundary  conditions 
are  easily  imposed,  and  because  it  reduces  to  a  two-dimensional  problem, 
thereby  reducing  the  computational  exercise. 
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We  continue  to  use  the   4>(x,t)   =  0  gauge.     The  known  time  and  longitudinal 
dependence  are  imposed  by  writing  the  vector  potential    and  current   induced  in 
the  walls  as 

A(x,t)   =  A(x±)   cos(o$t  -  kz)    =  A(l)   cos(wt   -  kz), 

(3.9) 


/(J)2  IT'' 

j    k    ~~    »  To   ~  ~' 


J(x,t)   =  J(i)  cos(ut  -  kz), 

c         a 

where  c  is  the  speed  of  light  in  the  waveguide  medium,  x  =  (x,y)  is  the  two- 
dimensional  transverse  position  vector,  and  where  A(x,t)  and  A(i)  are  real. 
Substitution  of  eq  (3.9)  in  eq  (3.2)  yields 

a+A     b+A     2  2 

S  =  C  J   dx   j   dy  ^  [A2(l)  +  A2(l)  ]  +  ^  A2(l) 
-A       -A  y 


(3.10) 


-  [(3y  Az(l))2+  (3x  Az(l))2+  (3x  Ay(l)  -  3y  Ax(l))2] 

+  2y  3(1)  •  A(l)}, 

C  =  4-  ]  dt  /  dz  cos2(ut  -  kz)  «  y~  j   dt  /  dz  sin2(wt  -  kz). 

In  writing  eq  (3.10)  we  have  assumed  that  the  material  in  the  waveguide  is 
isotropic  and  that  the  range  of  the  t  and/or  z  integration(s)  is  either  very 
long  or  an  even  number  of  cycles.  The  induced  current  term  and  the  limits  of 
the  transverse  integrations  require  explanation.  In  order  that  the  stationary 
point  of  S  yield  Maxwell's  equations,  the  variations  of  A  must  be  zero  on  the 
boundary,  which  in  our  two-dimensional  case  here  means  that  we  must  specify 
A  on  the  transverse  boundary.  In  order  to  be  able  to  specify  all  components 
of  A  we  choose  the  boundary  to  lie  a  few  skin  depths  within  the  conducting 
walls  where  the  field  can  safely  be  assumed  to  vanish.  Referring  to  figure  3- 
2,  the  x  integration  goes  from  x  =  -A  to  x  =  a  +  A,  and  the  y  integration  from 
-A  to  b  +  A,  where  A  =  N6,  some  suitable  number  of  skin  depths.  Then,  how- 
ever, the  currents  induced  in  the  walls  are  contained  within  the  integration 
volume.  Fortunately,  we  can  show  that  the  contribution  to  the  action  from  the 
volume  within  the  conductor  is  negligible,  and  we  can  write 
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S  -  C  J  dx  /  dy  \g   [Ay2(l)  +  A2(i)]  +^A|(l)  -  [(3  Az(l))2 

0       0 

2  2  (3"U) 

+  [3X  AZ(D)  +  ty  Ay(l)  -  3y  ax(i))  ]}  . 

We  still  need  to  impose  E.   =0  and  B     =  0  at  the  conductor  walls, 

tan        norm 

which  is  accomplished  by  the  requirements  that 

Az(0,y)  =  Az(x,0)  =  Az(a,y)  =  A2(x,b)  =  0, 

Ax(x,0)  =  Ax(x,b)  =  0,  (3.12) 

Ay(0,y)  =  Ay(a,y)  =  0. 

The  action  is  then  discretized  by  breaking  up  the  waveguide  into  a  rectangular 

grid  as  in  figure  3-3,  with  area  elements  centered  on  crosses  and  field  values 

defined  on  dots.  The  spacing  between  dots  is  Ax  =  a/N  ,  Ay  =  b/N  .  The  field 

x        y 

value  for  an  area  element  is  given  by  the  average  of  the  values  at  the  four 
corners  of  the  element.  For  derivatives,  the  average  of  th.e  two  appropriate 
differences  is  used,  e.g., 

d  A 


hnr)  H>  j4)  ■  2a7  tVi+1>  j+1>  -  V1'  ^ 


(3.13) 


+  Ax(1+l,j)   -  Ax(i,j)]. 


The  action  then  takes  the  form 


N         N 

S   -  C  Ax   Ay        j[  I        {(|-)2  DUa.e)   +  A(a-1,3)   +  A   (a,  3-1)   +  A   (a-1,3-1)] 

a=l      3=1  ™  y  y  y  y 

+   (|iT)2  [Ax(a'e)   +  Ax(a-1,3)   +  Ax(a,3-D   +  Ax(a-1,3-1)]2+   (^)2[Az(a,3) 
+  Az(a-1,3)   +  Az(a,3-D   +  Az(  a-1 , 3-1)]2  -    [j^f  [A^a.3)    -  Az(o,H) 


+  Az(a-1, 3)   -  Az(  a-1,3-1)  ]     -{-^)      [Az(a,3)   -  Az(a-1,3)  (3-14) 


+  A  fo,H)  -  A7(a-1,&-1)]- 


-  [-25J-  (Ay(a,3)  -  Ay(a-1,3)  +  Ay(a,3-1)  -  Ay(a-1,0-1)) 

-  J^   (Ax(a,3)  -  Ax(a,3-1)  +  Ax(*-1,3)  -  Ax(a-1,3-1)  )]2}   . 


We  next  wish  to  fix  the  field  values  at  a  small  number  of  measurement 
points  and  vary  all  the  fields  not  fixed  by  measurement  or  boundary  condition 
until  a  stationary  point  is  found.  Because  of  our  gauge  choice  and  the  fact 
that  we  are  only  considering  a  single  frequency,  a  measurement  of  E  is  a 
direct  measurement  of  A.   The  numerical  calculation  then  proceeds  as 
follows.  The  fields  A  are  fixed  at  measurement  points,  and  the  appropriate 
components  are  set  equal  to  zero  at  the  boundaries  (eq  (3.12)).  These  are  not 
allowed  to  vary  during  the  computation.  A  starting  field  configuration  is 
generated  according  to  the  prescription 


Ax(x,y) 


Ay(x,y) 


Az(x,y) 


B  (y)    m 

— I     A  (X  ,Y  ) 

x  n  n 


1 


R(x)  n=l 
B  (x)    m 


x  -  X 


B  (Y  )  ' 
x  n 


R(x)   n=l 


1         A  (X  ,Y  ) 
y  n  nJ 


x  -  X 


B  (X  )  ' 
y  rr 


N 


B  (y)  B  (x)    m 

— " I       A  (X  ,Y  ) 

R(x)      n=l   z  n  n 


1 


(3.15) 


B  (Y  )  B  (X  )  ' 
x  n   y  n 


where  R  is  just  a  normalization  factor  and  the  B  functions  enforce  the  bound- 
ary conditions, 

N 


R(x)  =  I 
n=l 


1 


x  -  X 


B  M  =  (i  .  e-2^/?)(l  -  e-2-3(^'/5)  , 


(3.16) 


B  (X)  =  (1  -  e-2-3x/5)(l  -  e-2-3(a-x)/e)  , 


17 


with  £  chosen  to  be  A/4.  The  vectors  X  =  (X  ,Y  1  are  the  measurement  points, 

n   v  n  n J  r  ' 

and  Nm  is  the  number  of  measurements.  The  initial  configuration  given  by  eqs 
(3.15)  and  (3.16)  is  an  arbitrary  choice,  but  it  does  have  the  desirable  prop- 
erties that  it  smoothly  interpolates  between  measurement  points,  obeys  the 
boundary  conditions,  and  is  not  the  correct  answer.  The  last  point  is  impor- 
tant because  the  question  is  whether  the  computation  can  find  the  correct  so- 
lution, not  whether  it  recognizes  the  answer  if  given  it. 

Having  generated  initial  values  for  all  the  A(i,j),  we  go  through  the 

grid  setting  the  values  of  A  at  each  point  (i,j)  equal  to  the  values  required 

to  make  3S/3A  (i,j)  =  0  given  the  current  values  of  A  at  neighboring  points, 
a 

Such  a  procedure  would  drive  us  to  a  nearby  extremum  if  only  there  were  one. 

In  either  the  continuous  (3.11)  or  discrete  (3.14)  form,  however,  the  action 

can  be  shown  to  have  no  (finite)  extremum.  The  numerical  procedure  will  then 

run  away,  given  enough  time.  (This  is  true  for  any  grid  size;  it  just  takes 
longer  for  finer  grids.) 

There  is  a  way  to  locate  the  stationary  point  nevertheless.  At  a  sta- 
tionary point,  S  should  change  very  little  as  we  pass  through  the  grid 
changing  the  A(i,j)'s,  and  so  we  plot  S  as  a  function  of  the  number  of  passes 
through  the  grid.  Figure  3-4  shows  the  result  for  a  =  0.8X,  b  =  0.4A.  The 
program  was  told  A  at  nine  points,  (Xn,Yn)  with  Xn  =  a/4,  a/2,  3a/4,  and  Yn  = 
b/4,  b/2,  3b/2.  A  40  x  20  grid  was  used.  The  slope  of  the  curve  is  quite 
small  from  N  -   40  to  N  ==  90;  the  approximate  stationary  point  is  somewhere  in 
this  range.  The  exact  point  doesn't  matter  much  since  the  fields  do  not 
change  much  in  this  range  of  N.  Figure  3-5  plots  the  fields  as  a  function  of 
x  for  a  few  values  of  y  for  N  =  65.  The  correct  solution  is  the  TE^q  mode, 
Ax  =  Az  =  °>  A  (x»y)  =  sin  ~  (solid  line  in  figure).  The  crosses  represent 
A  ,  the  circles  Ax.  The  computed  Az  is  zero  for  all  x  and  y  and  is  not 
plotted.  The  results  are  clearly  very  good  in  this  admittedly  simple  case, 
suggesting  that  the  method  holds  some  promise  for  practical  applications. 

We  have  also  performed  tests  in  which  the  number  and  positions  of  the 
measurement  points  were  varied.  Qualitatively,  the  results  are  about  what  one 
would  expect.  As  fewer  measurement  points  are  taken,  or  when  they  are  clumped 
together  or  all  far  from  the  center,  the  computed  solution  deteriorates.  A 
worst  case  is  only  one  measurement  point  (which  only  sets  the  scale  of  A  )  far 
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from  the  center.  For  waveguide  and  grid  as  above,  but  with  one  measurement 
point  at  (X,,Y,)  =  (0.1,0.1)  x,  the  action  behaves  as  in  figure  3-6.  Again  the 
general  location  of  the  approximate  stationary  point  is  apparent.  Taking  it 
to  be  at  N  =  55,  one  obtains  the  field  configurations  of  figure  3-7.  Although 
the  agreement  with  the  correct  solutions  is  not  so  good  as  it  was  with  nine 
measurements,  it  is  still  recognizable  and  would  be  quite  useful  if  we  only 
wanted  an  estimate  of  the  maximum  field  within  the  enclosure. 

This  simple  example  may  not  have  been  a  very  demanding  test,  but  the 
finite-element  action  approach  did  work;  and  it  worked  well  enough  to  bolster 
our  hope  in  the  eventual  practical  applicability  of  the  procedure. 

3.4  Problems  and  Prospects 

Although  the  waveguide  example  given  above  is  a  simple  problem,  the  suc- 
cess of  the  finite-element  action  approach  is  significant  nonetheless.  The 
program  did  find  a  solution  to  a  simple  two-dimensional  problem  without 
knowing  the  source(s).  Furthermore,  many  of  the  complexities  of  practical 
problems  are  little  more  than  technical  details  in  this  approach.  The  exten- 
sion to  three  spatial  dimensions,  for  example,  requires  more  computational 
time,  but  no  new  concepts  or  techniques.  Also,  because  the  symmetry  of  the 
waveguide  has  not  been  used,  extension  to  irregular  geometries  should  be  no 
more  difficult  than  specifying  the  geometry  for  the  computer  (provided  the 
geometry  is  not  such  that  it  requires  a  hopelessly  large  number  of  grid 
points).  Pieces  of  dielectric  material  can  be  included  simply  by  storing  the 
permittivity  at  each  grid  point,  e(i,j),  and  restoring  e(i,j)  to  eq  (3.15). 

The  biggest  unforeseen  difficulty  so  far,  that  the  stationary  point  is  a 
saddle  point,  has  been  overcome.  There  are  probably  more  elegant  and/or  effi- 
cient ways  of  finding  the  stationary  point,  but  we  have  demonstrated  that  at 
least  there  is_a_way  to  find  it.  It  is  possible  that  the  procedure  can  be 
modified  so  that  the  stationary  point  is  a  minimum,  but  we  do  not  know  how  at 
this  time,  and  it  is  not  necessary  as  long  as  saddle  points  can  be  located. 

All  this  is  not  to  say  that  we  are  home  free.  It  is  possible  that  in 
practical  problems  the  functional  structure  of  S  renders  the  stationary  point 
much  harder  to  locate.  In  addition,  there  are  a  number  of  technical  hurdles 
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to  be  overcome  before  this  technique  becomes  an  everyday  tool,  but  none  of 
them  is  obviously  insuperable.  In  a  practical  case  boundary  conditions  could 
cause  some  trouble.  Irregular  boundaries  are  no  problem  in  principle,  but 
they  will  entail  computational  complications.  A  more  serious  difficulty  is 
the  handling  of  open  boundary  conditions.  In  the  example  given,  the  volume 
was  completely  enclosed,  but  one  component  was  effectively  free  at  each 
wall.  When  all  three  components  are  free  over  some  part  of  the  boundary 
surface,  it  is  likely  that  some  measured  points  will  need  to  be  on  or  near 
that  surface  to  get  good  results.  In  the  waveguide  case  we  did  try  fixing  the 
value  of  all  components  of  A  on  the  boundary,  and  the  only  major  effect  was 
very  close  to  the  walls.  This  suggests  that  when  open  boundaries  are  present 
the  buffer  zone  (cf.  fig.  3-1),  wherein  we  solve  for  the  fields  but  do  not 
believe  the  results,  will  not  need  to  be  too  large. 

There  could  also  be  difficulties  associated  with  the  measurements  in 
practical  problems.  If  more  than  one  frequency  is  present,  then  a  measurement 
of  E  is  not  a  direct  measurement  of  A,  but  rather  of  3A/3t.   Consequently 
not  A  but  rather  its  derivative  must  be  held  fixed  in  the  computation. 
Similarly,  a  measurement  of  the  magnetic  field  would  fix  a  combination  of 
derivatives  of  A,  not  A  itself.  Another  difficulty  is  that  it  is  much  easier 
to  measure  the  peak  or  average  E  than  it  is  to  measure  E  at  one  time  t,  but 
the  latter  quantity  is  the  easier  to  include  in  the  calculation.  And,  of 
course,  eventually  we  need  to  confront  the  fact  that  the  measurements  are  not 
perfect,  and  uncertainties  must  be  included. 

Other  than  the  measurement  complications  just  mentioned,  complicated  time 
dependence  should  not  raise  major  difficulties.  It  increases  the  dimension  of 
the  problem—real  problems  will  be  four-dimensional  as  opposed  to  the  two- 
dimensional  example  we  did — but  it  is  handled  in  much  the  same  way  as  the 
spatial  dimensions.  It  increases  the  computing  time  required,  of  course,  but 
it  requires  no  new  developments. 

As  for  the  computation  itself,  no  great  effort  has  been  made  as  yet  to 
make  it  fast  or  efficient.  For  less  simple  geometries  a  different  (e.g., 
triangular)  gridding  system  would  be  more  versatile,  and  the  discretization  of 
the  integral  (3.11-14)  could  also  be  improved.  Such  refinements  and  sophisti- 
cations, however,  fall  in  the  fine-tuning  category.  The  more  immediate  task 
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is  to  advance  to  more  realistic  problems  and  try  to  handle  the  concomitant 
complications.     It  appears  that  this  action-based  approach  may  be  a  viable 
practical   method.     If  that  proves  to  be  the  case  we  would  have  a  very  powerful 
method  indeed;  but  considerable  work,  and  perhaps  even  some  ingenuity,  is  re- 
quired before  we  reach  that  point. 

4.     Scanning  Techniques 
4.1     Cylindrical   Scanning 

The  approach  that  is  presented  in  this  subsection  involves  cylindrical 
scanning.     The  ease  with  which  a  complex  EM  environment  can  be  scanned  by 
moving  a  highly  directional   probe  along  a  single  axis  for  different  azimuthal 
orientations  makes  cylindrical   scanning  very  attractive.     The  theoretical   de- 
velopment given  here  is  based  on  the  source-scattering  matrix  [13],  which  is 
very  analogous  to  the  well-established  plane-wave  matrix  approach  [14]. 

A  highly  directional   probe  antenna  and  its  coordinate  system  fixed  to  the 
probe  antenna  are  shown  schematically  in  figure  4-1.     Unknown  radiating 
sources  are  located  outside  a  cylindrical   volume,    p  >   p0.       The  EM  fields  in 
the  source-free  region   (p  <  p0)     will   be  expanded  in  a  complete  set  of  cylin- 
drical  eigenfunctions, 

00  00 

E(p,4>,z)   =       I         1      [b}(Y)  Mn   (p,M) 

n  ny 

f|=  — 00  —00 

+   b2   (y)    N      (p,<M)]  dY, 
n  nY 

H(p,<j>,z)   =  —  VxE  (4.1) 

0 

00  00 

mT~       I  1      My)   N      (P,<I>,Z)    +  b2M      (p,<J>,z)]  dY 

Zn  n  nY  n     nY 

U        D=_co        -oo 

I    =  /uTTT  «  377  a, 

0  0      0 


where  the  M  and  N     functions  are  given  by 


MnT(p,*,z)   =   Vx(Jn(KP)  e-jn<1)e-JYZ  ez) 
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=    p^J(Kp)  e     -  KJ'(Kp)   eJe"jnVJYZ 
Lp        n  p  n  <j)J 

NnY(p,*,z)  =  £  Vx  M^  (4.2) 

4  [-J*  J>P)  ep  -f=On(Kp)e;+  K2Jn(KP)ez]  e"H"^. 

All  quantities  with  a  hat  ( *)  over  them  denote  unit  vector;  eJa)  (w>0)  time 
dependence,  and  the  rationalized  mks  systems  are  used  throughout.  The  free- 
space  wave  number  k  is  defined  by  k  =  w/y0e0  =  "~  =  "T"  »  tne  axia^  Part  °f  tne 

C A 

wave  propagation  vector  is  denoted  y  and  K  =  A2  -  y2.  If  one  can  use  a  probe 
antenna  which  measures  the  transverse  (to  eJ  components  of  the  electric 
field,  then  the  b  (y)  can  be  determined  directly.  Namely,  using  the  ortho- 
gonality relationship  for  the  M  and  N  functions,  i.e. 

00  2  IT  „  «>  2  IT 

1       I     (I     xH.    ,)«e     dHzs0,     J       1     (Rnn  x  N.    ,  )   •  e     d*  dz  =  0 
J        n      ^  ny         n   y   '         p  '        h      K  nr         n  r    ;         p     Y 

-  ""  (4.3) 

L    So     [RnyxRn'Y^   -%d+dz=^Jn(KP)  J>p)   6_n>n,6(T+y). 
By  crossing  N  into  eq   (4.1),  one  gets  the  b*(y),   i.e. 


bn^>  '  4^"   (Kp0)  J1    (Kp0)     >       J       fl  x  E(p„.*.z)}   •  e*    d*  dz 

j         J     [-21  E  (p0,4>,z)  +  K*E  (Po,4,,z)]ejn*ejrZ  d*  dz. 


Similarly,  cross  M  into   (4.1)  to  yield  the  b2  (y), 

bn   ^  ■  4,^0   (KPo)     '  "/    £'.•♦•*>  *in^*  ^  ^ 

n       u     -oo     u 

Equations  (4.4,4.5)  give  the  modal  coefficients  in  terms  of  the  measured 
transverse  electric  field.  The  efficient  computation  of  the  double  integrals 
in  eqs  (4.4,4.5)  can  be  carried  out  by  use  of  the  sampling  theorem  and  a  fast 
Fourier  transform. 
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Consider  a  more  general  probe  antenna  whose  receiving  functions  R  (y)  are 
known  with  respect  to  a  cylindrical  coordinate  system  fixed  in  the  probe  and 
centered  on  itself.  The  output  of  the  probe  is  given  by 

b'o  =  I       I  J   R*  (Y)  a;s  (Y)  dY  ,  (4.6) 

s=l  n=-°°   -°° 

where  a'       are  the  modal   coefficients  of  the  J   (Kp)   modes  which  are  excited  by 
sources       existing       outside       the     probe     in     the       cylindrical       coordinate 
(  p'  ,  <f>'  ,z')   system  fixed  in  the  source. 


Specifically,  the  a     field  in  cylindrical   coordinates   (p'j^'jZ1)   is  given 


by 


00  00 


E'(p'.V.z')   -     I         /     [a;i(Y)  ^(p'.^'.z')   +  a'2(Y)   ^(p'.V.z'JldY. 

(4.7) 


00  00 


R'(p\*',z')  --±-     i      /    [b;MY)  ^(p'.f.z1)  +  b;i(Y)  ^(p'^'.z-jjdY 

"    o     n=_oo    —oo  '  ' 


Now  M„  is  defined  as  before, 
ny 

M^Ufcz)   =   VxJn(KP)e-jnVJTZ  ez 

=    [^nJn(KP)   -  u;   (KP)  e^]e-J'nV^ 

R^U.fcz)    =  17  V  x  M(1) 
ny  k  ny 


(4.8) 


4  U*  Jn^  e>^Jn^%+  K2jn^   eJe'^V^. 

Now,  regardless  of  whether  the  probe  antenna  or  EM  source  coordinate  system  is 
used,  the  EM  field  at  every  point  in  space  must  be  the  same. 

E'  Cp' »♦' »z<  )  =  E(p.4>»z) 

(4.9) 

H'tp'^'.z'  )  =  R(p,«|.,z] 
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when  (p',<t>,,zl)  and  (p,<)>,z)   refer  to  the  same  point  in  space.  But  as  the 
probe  antenna  scans  the  radiation  source,  the  source  coordinate  system  is 
merely  rotated  through  an  angle  <t>0,  about  the  z-axis,  and  translated  a  dis- 
tance z0  along  the  z-axis,  i,e. 


p'  =  p,  <j>'  =  <|>-<j>0  and  z'  =  z-z0  (4.10) 


and  eq  (4.9)  becomes 


E'  (p,  <t>-<J>0,  z-z0)  =  e(p><M) 

(4.11) 
H'  (p,<l)-  <|>0,  z-z0)  =  H(p,<j),z) 

i  i 

when  (E,H)  and  (E,H)  are  written  explicitly  in  terms  of  the  linearly  indepen- 
dent cylindrical  waves  in  eqs  (4.1)  and  (4.7).  The  orthogonality  relations  eq 
(4.3)  show  that  the  only  way  eq  (4.11)  can  be  satisfied  is  if 

*ln*(y)   =  b^(Y)e-jnVjnzo,   s-1,2.  (4.12) 

Thus  the  probe  receiving  equation  becomes 

•         00  00 

b'Uo.Zo)  =   I         I         J  R^(y)  b^Y)e-jn<t>0"jz°  dy  .       (4.13) 
s=l  n=-°<>  -» 

c 

Let  us  discuss  how  to  evaluate  some  characteristics  of  b  (y).  In  principle,  as 
the  probe  antenna  scans  the  EM  environment,  the  output  b'(<j>o,z0)  of  the  probe 
antenna  is  recorded  for 

0  <  <|>0  <  2tt  and  -~  <  z0  <  °°. 

That  is,  the  amplitude  and  phase  of  b'(<j>  ,z  )   are  measured  for  all   values  of 
<t>0  and  z  0.       In  practice,  the  zQ  scan  can  be  limited  to  some  finite  scan 
length,  since  b'0(  4>o»z  o)   ^s  assumed  negligible  outside  this  region.     Also,  data 
need  be  sampled  and  recorded  only  at  a  finite  number  of  measurement  points. 
Having  measured  b*  ( <t>0,z 0) »  the  Fourier  series  and  integral   of  eq   (4.13)  can  be 
immediately  inverted  to  yield  the  solution 
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I    Rs   (Y)  b*   (Y)  =   [Rj  (y)  bj  (y)  +  Rn2  (y)  b*  (Y)] 
„  .     n  n  n  n  n  n 

S=l 


2-rr 

1         r         r     ut,  .      _    n    Jn*n  JYZ( 


=  4^2     /       /    b'0U0,z0)  eJ"^eJ,to    d*0dz0  e  I  (y) 
4lT     -co       0  n(4.14) 

1  2 

Since  eq   (4.14)   involves  two  unknowns  b     (y)   and  b     (y),  two  linearly  inde- 
pendent scans  are  necessary  to  account  for  the  polarization  of  the  EM 
fields.     The  second  scan  produces  a  second  equation  to  complement  eq   (4.14), 
i.e., 

2       's  s 

I    RnS(Y)  b*   (y)  =  R^(y)  bi(Y)  +  Rn2(Y)  b2(Y) 

S=1  (4.15) 

1         °°      2X  jn<j)0     j'yzq 

=  — p     /       /      b'0  U0,z0)  e  e  d<j>0  dz0  =  I_'  (y) 

TV  _oo  U 

where  R  (y)  and  bo(<j>0,z0)  are  the  receiving  characteristic  and  the  output  of 
reoriented  probe,  respectively. 

Assuming  for  the  moment  that  the  receiving  function  (R1 1,  R'2)  of  the 
probe  antenna  is  known,  eqs  (4.14)  and  (4.15)  can  be  solved  immediately  for 
the  source  characteristic  b(Y),  i.e. 


bX(Y)  =  [Rn2(Y)  In(Y)  -  R2(Y)  In(Y)]/An(Y) 


b2(Y)  =  [RnX(Y)  i;(Y)  -  R^(Y)  In(Y)]/An(Y) 


(4.16) 


where 


An(Y)  =  R;2(Y)  RnMY)  -  R2(Y)  R'^Y)   .  (4.17) 


Equation  (4.16)  indicates  that  the  cylindrical  wave  expansion  coefficients 

b1  and  b2  for  the  EM  complex  field  environment  can  be  determined  provided  the 
n     n  r  r 

receiving  characteristics  of  the  probe  antenna  are  known.  There  are  several 
ways  to  determine  the  receiving  characteristics  of  a  general  probe  antenna  in 
terms  of  the  scalar  cylindrical  waves.  Since  it  is  beyond  the  scope  of  this 
paper,  it  will  be  discussed  in  a  future  paper.  But  it  is  very  easy  to  show 
that,  if  a  probe  antenna  is  an  ideal  electric  dipole,  i.e.,  a  probe  that  mea- 
sures the  transverse  components  of  the  electric  field,  then  eqs  (4.4)  and 
(4.5)  emerge  immediately  from  (4.16)  and  (4.17).  For  more  directional  probes, 
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a  different  expansion  may  be  more  useful,  and  we  next  turn  our  attention  to 
this  possibility. 

4.2  Spherical  (Directional)  Scanning 

If  one  has  a  probe  which  is  highly  directional,  this  feature  can  be  used 
for  directional  or  spherical  scanning  at  one  point,  a  suggestion  due  to  Chang 
and  Maley  [15].  The  directional  features  are  most  directly  exploited  by  a 
plane-wave  expansion.  Assuming  a  single  frequency,  the  electric  field  is 
written  as 

E(x,t)   =  E(x)   ejwt, 

3  -  -  (4.18) 

Zfo  =  /77T3e"Jk'X  l^   ■ 

(Ztt)3 

If  E  is  to  satisfy  Maxwell's  equations,  we  must  have  k2  =  <d2/c2,  allowing  us 
to  write   (for  large  volumes) 

E(Jc)  =  ek(k)  2w«(|-  k), 


It   is  necessary  at  this  point  to  set  forth  the  notational   conventions  for 
this  subsection  since  many  different  angles  will   be  encountered.     There  are 
two  coordinate  systems  of  interest,  one  fixed  with   respect  to  the  volume  of 
interest   (the  "lab"  coordinate  system),  and  one  fixed  with   respect  to  the 
probe  (the  "probe"  coordinate  system),  cf.  figure  4-2.     The  two  coordinate 
systems  have  a  common  origin,  the  location  of  the  probe   (assumed  to  be 
small).     The  z,  axis  of  the  probe  system  is  described  by  angles   e 0, <j> 0  in  the 
lab  system;  the   £  axis  of  the  probe  system  is  chosen  to  lie  in  the  z  x   z, 
plane.     In  the  course  of  the  scanning   e0  and   $$  vary,  of  course.     The  direc- 
tion of  incidence  of  a  plane  wave  with   respect  to  the  probe  system  will   be 
denoted   6',<j>'.  This  same  direction  as  seen  in  the  lab  system  is  called 
8,<j).  The  angles    0,<j>  depend  not  only  on   6',$'   but  also  on   9q»<I>q.     And  con- 
versely  e'   =  e'  (8,<t>,e0,<j>0),   y   =  <j>'  ( e, cj>, 60, 4>o)  •  A  straightforward  exercise 
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with  rotation  matrices  yields  the  relations 

cos  0'  =  sin  e0  sin  8  cos(  (|>- 4>0)  +  cos  8o  cos  e» 

(4.20) 
cot  <j>'  =  cot(  <j»-<t>o)  cos  e0  -  sin  e0  cot  0  csc(<t>-<t>0). 

For  the  sake  of  simplicity  and  clarity  in  this  presentation,  we  limit 
ourselves  to  the  case  of  a  scalar  field.  The  response  of  the  probe  to  a 
single  plane  wave  incident  at  angles  0',<j>'  in  the  probe  coordinate  system  will 
be  the  product  of  the  probe  acceptance  at  that  angle  (A),  the  amplitude  per 
solid  angle  of  the  wave  (eM,  and  the  solid  angle, 

dR  =  A(e'  ,V  )  ej,(0',<t>')  da'.  (4.21) 

But  the  plane  wave  component  from  0' ,<j>'in  the  probe  system  is  the  same  as  the 
component  from  0(  0' ,<(>',  0o»<j>o)  >  <K  e' »  ♦' »  9o»  H)  in  tne  ^aD  system, 

e'k(0\<j>')  =  ek(0(0' ,<t»' ,0o,<t>o),  4>(e' ,4>' ,e0,<|>o)).        (4.22) 

Using  this  fact  and  integrating  the  differential  response  (4.21)  over  all  in- 
cident angles  in  the  probe  system,  we  obtain  for  the  total  response 


1  TT 

R(  ©o»  4>0)   =     /     d  cos    0'      /       d<j>'   A(  8' ,  4)' )  (4.23) 

-1  -TT 

e^  (0(  0' ,  <j>' ,  0o>  <J>o)  >   4»(  e'  s^1  »©o»4>o)  J- 

For  purposes  of  numerical    inversion,  we  find  it  more  convenient  to  rewrite  eq 
(4.23)   as  an  integration  over  lab-system  angles, 

1  TV 

R(e0»*o)  =  /     d  cos  0  /  dcf.  A(e'  ( e,  <j>,  e0, <t>0) ,  <t>'  ( e, <j>, e0, <j>0) )  ek(  e,<t>).   (4.24) 

-1  -TT 

If  the  probe  is  perfectly  directional  A  is  just  a  product  of  delta  functions 
of  0  -  0O  and  4)  -  <j>0;  and  the  response  at  0o,<j>o  measures  e.  ( e0, <j>0)  direct- 
ly.    For  imperfect  directional   probes,  a  little  more  work  is  required.     The 

probe  response  is  recorded  for  some  number  N  of  scan  angles 

r  r  meas 
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9,-»*:»i   =  ^Nmeas*     We  must  invert   (4.24)  to  obtain  e.  (6,9),  which  then  deter- 
mines the  field  throughout  the  volume  from  eq   (4.19). 

To  invert  eq   (4.24)  we  use  a  pulse  expansion  of  ek   (cos0,<j>)  and  point 
matching.     Assume  for  simplicity  that  the  scanning  was  done  at  regular  inter- 
vals   6n  in  cose  and   6    in   <t>.     Let 

o  <p 

_      Ne    \ 
e,(e,<j>)  =       I      I      enm  n(cos   e,<t>), 


'k       nm 
n=l  m=l 


Ne=2/6e>  %=2^V 


(4.25) 


where    n    (cos    e,<b)   is  the  unit  pulse  function  centered  at 
nm  ,  r, 

cos    e  =  -1   +(n  -  2-)6Q,    <j>  =   (m  -  j)6 


n     (cose, 9)  =1        -1  +  (n  -  1)6Q  <  cose  <  -1  +  n6, 


and   (m-l)6A  <   <t>  <  m5x,  (4.26) 

9  9 

=  0  otherwise. 


For  each  of  the  Nmeas  measurement  directions  we  get  an  equation  of  the  form 

Ne      Nd> 

R(V»i)  -    I      I     ekm  An>i>  M> 

n=l       m-1  (4#2?) 

-l+n6  mS 

Ann,(9,f«.  )=       /    u    d  cose      J*    d<t>  A  (e'(e,  ^e. ,  $ .),  <>•  ( 0,  <j>.e ,,*_•)). 
1      '       -l+(n-l)6e  (m-1) 6^  n      ' 

Because  the  measurement  angles  coincide  with  the  centers  of  the  pulses,  it  is 

a  simple  matter  to  relabel  indices  and  force  eq  (4.27)  into  the  form  of  a 

simple  matrix  equation.  The  m  and  n  indices  are  combined  into  one  index 

running  from  1  to  N     =  NnxN,, 
3         meas    e    9' 

N  +1 
ai  =  eii    a2  =  ei2        a  9   =  p21 
9k  "  ek  '     9k  "  ek  9k     ek  ,  ...  , 

a\   =  ^l(Qv^),af  =  A12^.,*.),  ...  ,  (4.28) 
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6-   =  R(e.  ,  $. ). 
Then  eq   (4.27)  takes  the  form 


N 
meas 


*_1  (4.29) 

gJk  =   I  (a1)0'1    r 

The  e.  's  are  then  substituted  into  eq  (4.25)  to  yield  the  approximation  for 
ek(e,<t>),  which  in  turn  is  substituted  into  eq  (4.19)  to  yield  the  field  E(x). 

We  thus  can  obtain  an  estimation  of  E(x)  even  with  only  a  few  scanning 
angles,  though  of  course  the  approximation  becomes  better  as  more  measurements 
are  made.  It  remains  to  be  seen  how  large  the  angular  increments  in  the  scan- 
ning can  be  in  practical  applications.  The  expansion  of  eq  (4.19)  will  only 
be  valid  for  |x|  less  than  the  distance  to  the  first  source  or  scatterer,  but 
the  method  still  offers  the  promise  of  an  estimate  of  the  field  throughout  a 
volume,  from  a  limited  number  of  measurements  at  just  one  point. 

5.  Summary 

As  more  sources  contribute  to  ambient  electromagnetic  fields,  and  as  the 
electronic  devices  operating  in  these  ambient  fields  become  more  numerous, 
sensitive,  and  important,  the  problem  of  efficiently  measuring  and  charac- 
terizing complicated  electromagnetic  environments  is  becoming  increasingly 
acute.  We  have  discussed  three  different  approaches  to  the  problem,  outlining 
the  foundation,  present  status,  and  direction  of  future  development  of  each. 
Because  there  is  so  little  previous  work  on  characterization  of  complicated 
electromagnetic  environments,  any  new  method  tends  to  require  completely  new 
development  beginning  from  the  basics,  and  as  a  result  progress  can  be  rather 
slow. 

Each  method  discussed  is  quite  promising  in  the  sense  that  each  appears 
to  have  a  reasonable  chance  of  actually  working,  and  each  will  be  very  useful 
if  it  does  work.  However,  each  method  also  requires  further  development  and 
work.  That  work  is  in  progress. 
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Appendix 

In  order  to  locate  a  stationary  point  of  the  action  we  need  to  know 
whether  we  are  looking  for  a  maximum,  a  minimum,  or  neither.  In  addition,  if 
it  is  an  extremum,  we  would  like  to  know  whether  it  is  a  global  or  only  a 
local  extremum.  It  is  often  assumed  that  it  is  a  minimum,  but  our  experience 
indicates  that  exceptions  are  neither  so  rare  nor  so  pathological  as  one  might 
expect. 
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For  mechanical  systems  the  problem  is  addressed  in  Whittaker's  classic 
treatise  [9],  (Note  that  Whittaker  uses  the  word  "action"  for  a  different 
quantity  than  we  do.   It  is,  however,  a  difference  of  nomenclature  not  sub- 
stance.) The  type  of  stationary  point  which  S  has  depends  on  the  "kinetic 
focus"  of  the  original  point,  defined  as  follows.  Consider  a  point  A  on  an 
actual  trajectory;  and  let  another  actual  trajectory  pass  through  A,  at  some 
small  angle  with  respect  to  the  first.  If  the  two  trajectories  intersect 
again,  say  at  point  B,  then  the  kinetic  focus  of  A  is  the  limiting  value  of  B 
as  the  angle  between  the  two  trajectories  at  A  goes  to  zero.  The  relevant 
result  then  is  that  the  stationary  point  of  the  action  S  is  a  minimum  if  the 
final  point  in  the  integration  occurs  before  the  kinetic  focus  of  the  initial 
point;  if  it  occurs  after  the  kinetic  focus  the  stationary  point  is  neither  a 
minimum  nor  maximum  but  rather  a  saddle  point. 

To  make  sense  of  the  preceding  verbiage,  an  example  is  in  order. 
Whittaker  considers  a  particle  on  a  sphere  with  no  forces  acting.  For 
variety,  and  because  it  is  a  little  more  relevant  for  us,  we  consider  a  simple 
harmonic  oscillator  such  as  a  mass  (m)  on  a  spring  (spring  constant  k).  The 
action  is 

S  =f  J   dt  [x(t)2-£x(t)2],  (A.l) 

ti 
where  x(t)  is  the  position  at  time  t.  The  kinetic  focus  for  a  given 
x.  =  x(t.)  is  easily  determined  from  our  knowledge  of  the  solutions.  Since 
solutions  of  different  amplitude  all  have  the  same  period,  if  two  trajectories 
(solutions)  are  both  at  x^  at  time  t^,  then  one  half  period  later  both  will  be 
at  x  =  -x..  Consequently,  the  kinetic  focus  of  t^  occurs  at 
t  =  t.  +  T/2  =  t.  +  ir/rn/k  .  We  therefore  expect  the  stationary  point  of  S  will 
be  a  minimum  for  tf  <  t.  +  T/2  and  a  saddle  point  if  t-  >  t.  +  T/2. 

We  can  verify  this  expectation  by  direct  calculation.   If  we  vary  x(t) 
and  set  the  first  variation  of  S  equal  to  zero  we  just  obtain  the  equation  of 
motion.  Considering  the  second  order  variation  of  S,  we  have 

§<JS-  j   dt  [<*)*-*<*)*] 

*1  (A.2) 
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tf 


=    [6x    6x]^f  -      /       dt    6x    [Sx  +  -  Sx], 


,      t. 

The  first  term  vanishes  because   6x(t)   =  0  at  the  end  points  t  =  t.,  tf.  This 
fact  also  allows  us  to  expand   6x(t)   in  a  discrete  Fourier  sine  series, 

00 

6x(t)   =        I         Sx     sin    [n    ir(t-t.)/At]   , 

n_1  (A.3) 

At    =  tf    -  t  .    . 

Substituting  eq   (A.3)   in  eq   (A. 2)  yields 

I     SY2S  =#"       I         [(-Tt)2"-]   Sx2     .  (A.4) 

m       x         2         Ln         lKt&J         mJ       n  v        ' 

n=l 

Since  the     original    6x(t)  was  arbitrary,  so  too  are  all   the   Sxn's.     Therefore, 
62S  will   always  be  positive  only  if    [(mr/At)2  -  k/m]  >  0  for  all   n.     This  re- 
quires 

At   <  *  ^  =  IT*  (A*5) 

For  At  >  T/2,  62S  can  be  positive  or  negative  depending  on  the  choice  of 
6x(t).  Therefore,  as  expected,  the  stationary  point  is  a  minimum  if  tf  is 
earlier  than  the  kinetic  focus,  and  a  saddle  point  if  tf  is  after  the  kinetic 
focus.  (If  tf   occurs  at_  the  kinetic  focus,  62S  >  0.) 

Transferring  these  results  to  the  electromagnetic  problems  of  interest 
requires  further  work;  but  there  are  reasons,  both  heuristic  and  empirical, 
for  expecting  the  stationary  point  to^be  a  saddle  point.  One  heuristic  reason 
is  that  we  are  formally  considering  J  dt  L,  so  that  tf  must  be  later  than  the 
kinetic  focus  provided  it  exists.  Secondly,  in  the  mechanical  case  the  tran- 
sition from  minimum  to  saddle  point  is  related  to  the  nonuniqueness  of  the 
solution.  For  the  harmonic  oscillator  example,  for  tf  <  t.  +  T/2  specifying 
t.,  tf,  x(t.)»  and  x(tf)  specifies  a  unique  solution  (i.e.  amplitude  and 
phase).  However,  for  tf  =  t^  +  T/2  there  is  a  continuum  of  solutions  possible 
for  a  given  t-j,  tf,  x(t.)»  x(tf).  This  is  no  coincidence  since  the  kinetic 
focus  occurs  at  the  (limit  of  the)  intersection  of  real  trajectories  passing 
through  the  initial  point.  It  is  then  reasonable  to  infer  that  the  stationary 
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point  ceases  to  be  a  minimum  when  tf  exceeds  the  first  time  at  which  it  is 
possible  to  have  more  than  one  solution.  The  nonuniqueness  of  the  solutions 
in  the  electromagnetic  problems  we  are  attacking  was  noted  earlier,  and  so  we 
have  another  qualitative  reason  to  expect  a  saddle  point. 

Finally,  there  is  the  empirical  evidence:  in  the  example  problem  pre- 
sented in  subsection  (3.3),  the  stationary  point  is  shown  to  be  a  saddle 
point.  This  is  a  potential  disaster  computationally,  but  it  will  prove  pos- 
sible to  locate  the  saddle  point. 
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Figure  2-1.       Cumulation  amplitude  probability  distribution  of  EM  noise  in  a 
mine. 
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Figure  2-2.       Average  crossing  rate  of  EM  noise  in  a  mine, 
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Figure  2-3.       Cumulative  amplitude  probability  distribution  of  the  building's 
attenuations.^-'     The  straight  line  represents  a  log-normal   distribution. 
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Figure  2-4.   Cumulative  amplitude  probability  distribution  of  fields  in  a 
rectangular  cavity.'-7^  The  straight  line  represents  a  lognormal  distribu- 


tion. 
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Figure  3-1.   Schematic  division  of  total  volume  into  volume  of  interest  (Vj), 
buffer  volume  (VB),  and  remaining  volume  containing  all  sources  (V$). 
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Figure  3-2.   Rectangular  waveguide  dimensions  and  axes.  Also  depicted  is  the 
boundary  of  the  volume  over  which  the  action  integral  extends. 
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Figure  3-3.       Lattice  of  points  for  finite-elment  action, 
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Figure  3-4.   Action  as  a  function  of  N,  the  number  of  passes  through  the  grid 
in  the  numerical  computation,  for  nine  measurement  points. 
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Figure  3-5.       Fields  Ax  and  A    as  functions  of  x  at  fixed  y,  after  65  passes 
through  the  grid  in  the  numerical   computation  with  nine  measurement 
points.     Solid  line  is  correct  result  for  Ay.     Correct  result  for  Ax  is 
zero,     a)  y  =  0.2A    b)  y  =  O.lA    c)  y  =  0.05A 
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Figure  3-6.       Action  as   function  of   N,  for  one  measurement   point. 
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Figure  3-7.       Computed  fields  as  functions  of  x  at  fixed  y,  after  55  passes 
through  the  grid,  with  one  measurement  point.     Solid  line  is  correct  result 
for  Au.     Correct   result   for  A     is   zero,     a)  y  =  0.2X     b)   y  =  0.1X 
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c)  y  =  0.05X 
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Figure  4-1.   Coordinate  system  for  probe  antenna 
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Figure  4-2, 
ni  ng. 


Coordinate-system  conventions   for  spherical    directional    scan- 
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